Symbiodiniaceae community dynamics in Pocillopora colonies from Uva Island reef during 2014-2016
Ana Palacio-Castro
March 29, 2020
- 1. Symbiont to host (S/H) cell ratio and Durusdinium proportion
- 2. qPCR Data Cleaning and quality control
- 3. Field/sample information
- 4. Community classification by genera presence/dominance
- Initial qPCR BASED analysis (No ITS2 or ORF information included)… skip to 5b?
- Symbiodiiaceae ITS2 type + ORF lineage data analysis
- Added Pocillopora ORF lineage and Symbiodiniaceae ITS2 type hosted in 2014
- 5b. Changes in the number of colonies dominated by each genus
- 6b. Changes in colony proportions based on ITS2 profiles (2014) and genus dominance
- 7b. Changes in total Symbiodiniaceae abundance (bleaching): qPCR + ORF
- 8b. Changes in Cladocopium abundance
- 9b. Changes in Durusdinium abundance
- 10. Supplmenmtary C/H D/H dynamics
- 5c. Changes in the number of colonies dominated by each genus considering mortality
This file analyzes the changes in S/H cell ratio (bleaching and symboint shuffling) in Pocillopora spp.: * 29 colonies were sampled in 2014 (Before ENSO), 2015 (during a 1st peak of heat) and 2016 (after a 2nd peak of heat) * These samples were collected at Uva Island, Gulf of ChiriquÃ, Panama.
# Libraries
library(plyr)
library(tidyverse) # install.packages('tidyverse')
# library(devtools) # install.packages("devtools")
# devtools::install_github("jrcunning/steponeR")
library(steponeR)
library(scales)
library(gridExtra)
library(lme4)
library(lmerTest)
library(multcomp)
library(emmeans)
# Default ggplot settings
ggthe_bw<-theme(plot.background=element_blank(),
panel.grid.major.y = element_blank(),
panel.grid.major.x = element_blank(),
panel.grid.minor.x = element_blank(),
panel.grid.minor.y = element_blank(),
axis.line = element_line(colour = "black"),
legend.box.background = element_rect(),
panel.background =element_rect(fill = NA, color = "black")
)+
theme_bw()1. Symbiont to host (S/H) cell ratio and Durusdinium proportion
Get the raw data for all Pocillopora samples applying R.Cunning steponeR function
- Get list of plate files to read data
- Calculate the ratios
- Export results
Pdam.plates <- list.files(path="data", pattern=".csv", full.names=T)
Pdam.plates## [1] "data/Chi-Pdam_Plate2.csv" "data/Chi-Pdam_Plate3.csv"
## [3] "data/Chi-Pdam_Plate4.csv" "data/Chi-Pdam_Plate5.csv"
## [5] "data/Chi-Pdam_Plate6.csv" "data/Gal-Pdam_Plate13.csv"
## [7] "data/Gal-Pdam_Plate14.csv" "data/Gal-Pdam_Plate15.csv"
## [9] "data/Gal-Pdam_Plate16.csv" "data/Gal1_Dprobes_data.csv"
## [11] "data/Pan-Pdam_Plate10.csv" "data/Pan-Pdam_Plate7.csv"
## [13] "data/Pan-Pdam_Plate8.csv" "data/Pan-Pdam_Plate9.csv"
## [15] "data/Pdam_Plate11-Deep.csv"Pdam.Out <- steponeR(files=Pdam.plates, target.ratios=c("C.Pdam", "D.Pdam"),
fluor.norm=list(C=3.798, D=0, Pdam=5.416),
copy.number=list(C=9, D=1, Pdam=3),
ploidy=list(C=1, D=1, Pdam=2),
extract=list(C=0.813, D=0.813, Pdam=0.982))
# Target ratio results
Pdam<-Pdam.Out$resultGet colony and year information from sample labels
# Parse sample names, times, treatments, colonies, plates
sample.year <- rbind.fill(lapply(strsplit(as.character(Pdam$Sample.Name), split="_"),
function(X) data.frame(t(X))))
colnames(sample.year) <- c("Spp", "Colony","Year")
Pdam <- cbind(sample.year, Pdam[,-1])
# Keep unique sample ID to check reruns
Pdam$Sample<-paste(Pdam$Year,Pdam$Colony, sep="_")
Pdam$ID<-paste(Pdam$Sample,Pdam$File.Name, sep=".")Calculate total S/H ratio and D/C ratio and Proportion of Durusdinium in each sample
# 0. If Clade only detected in one technical replicate, set its ratio to NA (becomes zero)
Pdam$C.Pdam[which(Pdam$C.reps==1)] <- NA
Pdam$D.Pdam[which(Pdam$D.reps==1)] <- NA
# 1. Rename cols and make NA=0
colnames(Pdam)[which(colnames(Pdam) %in% c("C.Pdam", "D.Pdam" ))] <- c("C.SH", "D.SH")
Pdam$C.SH[is.na(Pdam$C.SH)] <- 0
Pdam$D.SH[is.na(Pdam$D.SH)] <- 0
Pdam <- Pdam[, c("Colony", "Year", "C.SH", "D.SH","Sample", "ID")]
# 2. Get the ratios and log 10
# Total ratio
Pdam$tot.SH <- Pdam$C.SH + Pdam$D.SH
Pdam$logTot.SH <- log10(Pdam$tot.SH )
# C ratio
Pdam$logC.SH <- log10(Pdam$C.SH)
# D ratio
Pdam$logD.SH <- log10(Pdam$D.SH)
Pdam$logTot.SH[which(Pdam$tot.SH==0)] <- NA
Pdam$logC.SH[which(Pdam$C.SH==0)] <- NA
Pdam$logD.SH[which(Pdam$D.SH==0)] <- NA
# Clade Proportion
# D Proportion
Pdam$D.Prp<-(Pdam$D.SH/Pdam$tot.SH)
# C Proportion
#Pdam$C.Prp<-(Pdam$C.SH/Pdam$tot.SH)
hist(Pdam$D.Prp)
#Dproportion<-Pdam[(Pdam$D.Prp<0.99 & Pdam$D.Prp>0.01),]
#hist(Dproportion$D.Prp)Note: Most of the colonies are dominated by one symbiont genus. Even proportions of multiple symbionts are rare.
2. qPCR Data Cleaning and quality control
## 1. Check and remove NTC wells
ntc <- Pdam[which(Pdam$Colony=="NTC"), ]
Pdam <- droplevels(Pdam[!rownames(Pdam) %in% rownames(ntc), ])
## 2. Check and remove NoHSratio samples
NoHSratio <- Pdam[which(Pdam$tot.SH==0), ]
Pdam <- droplevels(Pdam[!rownames(Pdam) %in% rownames(NoHSratio), ])
## 3. Chose bw samples ran more than once
ReRunA <- Pdam[duplicated(Pdam$Sample),]
n_RunA <- data.frame(table(Pdam$Sample))
colnames(n_RunA)<-c("Sample","RanA")
Pdam<-join(Pdam, n_RunA, type = "left")
DuplicatesA <- Pdam[(Pdam$RanA>1),]
# Remove bad replicates
ToRem1<-read.csv("ToRemove10-28-16.csv") # 11/1/16
Pdam<-Pdam[!(Pdam$ID %in% ToRem1$ID),]
n_RunB <- data.frame(table(Pdam$Sample))
ReRunB <- Pdam[duplicated(Pdam$Sample),]
colnames(n_RunB)<-c("Sample","RanB")
Pdam<-join(Pdam, n_RunB, type = "left")
# List of dupplicated samples, should have 0 rows now
DuplicatesB <- Pdam[(Pdam$RanB>1),]
# 4. Check and remove samples with high SD or late amplification
StDe1.5 <- Pdam[which(Pdam$Pdam.CT.sd>1.5), ]
Pdam<-Pdam[!(Pdam$ID %in% StDe1.5$ID),]
# End of qPCR data cleaning3. Field/sample information
Get collection information and select only samples from Uva Island that were sampled in 2014, 2015 and 2016
# 1. Labeling and Factors
## Import treatment info from Field file and mortality / ITS2 / ORF data
Info<-read.csv("SampleInfo.csv", header=TRUE)
Mortality<-read.csv("MortalityData.csv", header=TRUE)
Mortality$ITS_2<-factor(Mortality$ITS_2, level=c("C1d", "C1b-c", "D1"))
## Merge SH ratios and sample info
Pdam<-join(Pdam, Info, by = "Sample",
type = "left", match = "all")
#summary(Pdam)
# 2. Uva Data
# Write how many years a colony was sampled
Years <- data.frame(table(Pdam$Colony))
colnames(Years)<-c("Colony","Times")
Pdam<-join(Pdam, Years, type = "left")
# Select data only from Uva Island
data.Uva<-Pdam[which(Pdam$Location=="Uva Island"), ]
# Select colonies that were sampled three times (2014, 2015 and 2016)
data.Uva<-data.Uva[which(data.Uva$Times>2),]
# Remove "Ross" samples (colony in deep location (~50feet), not Uva reef)
# and 259 as well as 264 that were mislabeled in at least one of the sampling points
data.Uva<-filter(data.Uva, Colony!="Ross")
data.Uva<-filter(data.Uva, Colony!="259-1")
data.Uva<-filter(data.Uva, Colony!="264-1")
data.Uva<-data.Uva[,c("Colony","Year.2","C.SH","D.SH","tot.SH",
"logTot.SH","D.Prp", "logC.SH", "logD.SH",
"Year","Sample", "Reef", "Spp", "Depht")]
data.Uva$Year<-factor(data.Uva$Year, levels=c("14", "15", "16"))
data.Uva$Colony<-factor(data.Uva$Colony)
#summary(data.Uva)
#3. Merge SH ratios and moratlity info
data.Uva<-join(data.Uva, Mortality, by = "Sample",
type = "left", match = "all")
summary(data.Uva)## Colony Year.2 C.SH D.SH
## 256-1 : 3 Min. :2014 Min. :0.000000 Min. :0.00000
## 257-1 : 3 1st Qu.:2014 1st Qu.:0.000000 1st Qu.:0.00000
## 260-1 : 3 Median :2015 Median :0.003005 Median :0.01442
## 262-1 : 3 Mean :2015 Mean :0.019174 Mean :0.02071
## 278-1 : 3 3rd Qu.:2016 3rd Qu.:0.011841 3rd Qu.:0.03001
## 284-1 : 3 Max. :2016 Max. :0.690021 Max. :0.18875
## (Other):69
## tot.SH logTot.SH D.Prp logC.SH
## Min. :0.00146 Min. :-2.8356 Min. :0.0000 Min. :-5.0627
## 1st Qu.:0.01294 1st Qu.:-1.8888 1st Qu.:0.0000 1st Qu.:-2.5001
## Median :0.02266 Median :-1.6448 Median :0.8620 Median :-2.0971
## Mean :0.03988 Mean :-1.6539 Mean :0.5411 Mean :-2.1808
## 3rd Qu.:0.03579 3rd Qu.:-1.4463 3rd Qu.:1.0000 3rd Qu.:-1.7015
## Max. :0.69002 Max. :-0.1611 Max. :1.0000 Max. :-0.1611
## NA's :32
## logD.SH Year Sample Reef
## Min. :-5.7262 14:29 Length:87 Coral Collection Point:36
## 1st Qu.:-1.8381 15:29 Class :character Shallow Millepora :18
## Median :-1.6245 16:29 Mode :character M.intrincata Shallow :12
## Mean :-2.0059 Ridley Rock Shallow : 9
## 3rd Qu.:-1.4640 Gardineroseris city : 6
## Max. :-0.7241 Bw Gard City and 4x5 : 3
## NA's :29 (Other) : 3
## Spp Depht Mortality ID_Glynn ID_Palacio
## P.dam:51 Min. : 9.00 Almost total: 7 P. cap:12 P. dam:24
## P.eff: 0 1st Qu.:15.00 No :71 P. dam:27 P. eyd: 6
## P.ele:36 Median :17.00 Partial : 9 P. ver:48 P. ver:57
## P.eyd: 0 Mean :17.79
## P.mea: 0 3rd Qu.:20.00
## P.ver: 0 Max. :32.00
##
## ORF ITS_2
## Type_1:69 C1d :18
## Type_3:18 C1b-c:20
## D1 :49
##
##
##
## data.Uva$Mortality<-factor(data.Uva$Mortality, levels =
c("No", "Partial", "Almost total"))ORF types versus morphological ID
# Summary ID *Pocillopora* by Glynn
Summary_ID_Glynn<-data.Uva %>%
group_by(Year.2, ORF, ID_Glynn) %>%
dplyr::summarise(total.count=n())
Summary_ID_Glynn # Summary ID *Pocillopora* by Palacio
Summary_ID_Glynn<-data.Uva %>%
group_by(Year.2, ORF, ID_Palacio) %>%
dplyr::summarise(total.count=n())
Summary_ID_Glynn 4. Community classification by genera presence/dominance
Symbiodiniaceae genera detected in the samples (Only-C, Only-D, Mixed communities)
# Add a column for the genera detected in each sample using qPCR
data.Uva$Community<-NULL
data.Uva$Community[which(data.Uva$D.Prp>=0)] <- "Mixed"
data.Uva$Community[which(data.Uva$D.Prp==0)] <- "Cladocopium-only"
data.Uva$Community[which(data.Uva$D.Prp==1)] <- "Durusdinium-only"
data.Uva$Community<-factor(as.character(data.Uva$Community),
levels=c("Cladocopium-only", "Mixed", "Durusdinium-only"))
# Summary genera *Pocillopora*
Summary_Community<-data.Uva %>%
group_by(Year.2, Community) %>%
dplyr::summarise(total.count=n())
Summary_Community # Summary genera detected by ORF type
Summary_Community<-data.Uva %>%
group_by(ORF, Year.2, Community) %>%
dplyr::summarise(total.count=n())
Summary_CommunityDominant Symbiodiniaceae genus in each sample (C-dominated, D-dominated)
In general, Pocillopora colonies are dominated by one Symbiodiniaceae type that accounts for more than 80-90 of the symbiont community. However, sometimes intermediate proportions are common when the symbiont community is transitioning from being dominated by one algal type to other (e.g. during heat stress).
For this data set, 2 of the colonies sampled in 2014 had relativelly even proportions of Cladocopium_C1bc and Durusdinium glynnii (Colony 267_2014, D proportion = 0.47; Colony 269_2014, D proportion = 0.50). This made difficult their classification into one domminant type for statistical testing and further plotting. Since these 2 colonies were later dominated by Durusdinim glynnii in 2015, we decided to include them in the D-dominated classification, eventhoug in 2014 there were not such a dominance (See figure 2a).
# Add a column for the dominant genus
data.Uva$DominantCommunity<-NULL
data.Uva$DominantCommunity[which(data.Uva$D.Prp>=0.45)] <- "Durusdinium-dominated"
data.Uva$DominantCommunity[which(data.Uva$D.Prp<0.45)] <- "Cladocopium-dominated"
data.Uva$DominantCommunity<-factor(as.character(data.Uva$DominantCommunity),
levels=c("Cladocopium-dominated","Durusdinium-dominated"))
# Summary dominant genera *Pocillopora*
Summary_Community2<-data.Uva %>%
group_by(Year.2, DominantCommunity) %>%
dplyr::summarise(total.count=n())
Summary_Community2 # Summary dominant genera by ORF type
Summary_Community3<-data.Uva %>%
group_by(ORF, Year.2, DominantCommunity) %>%
dplyr::summarise(total.count=n())
Summary_Community3Set initial (2014) and final (2016) communities for each colony
data.Uva$ITS_2<-factor(data.Uva$ITS_2, level=c("C1b-c", "C1d", "D1"))
# Subset Aug 2014, Aug 2015 and April 2016 samples
D2014.data<-subset(data.Uva, Year.2==2014)
D2014.data<-D2014.data[,c("Colony","D.Prp","Community","DominantCommunity", "ITS_2")]
colnames(D2014.data)<-(c("Colony","D14.D.Prp","D14Community", "D14Dominant", "ITS2_14"))
data.Uva<-left_join(data.Uva, D2014.data, by=c("Colony")) # Add initial community to all samples
D2016.data<-subset(data.Uva, Year.2==2016) # D1a under control temperature
D2016.data<-D2016.data[, c("Colony","D.Prp","Community","DominantCommunity", "ITS_2")]
colnames(D2016.data)<-(c("Colony","D16.D.Prp","D16Community", "D16Dominant", "ITS2_16"))
data.Uva<-left_join(data.Uva, D2016.data, by=c("Colony")) # Add final community to all samplesInitial qPCR BASED analysis (No ITS2 or ORF information included)… skip to 5b?
5a. Changes in the number of Pocillopora colonies dominated by each genus
To test the statistical significance of changes in the number of colonies dominated by Durusdinium versus Cladocopium we considered inclussion and exclussion of these 2 initially even colonies, and to place them in the D-dominated or C-dominated categories in 2014.
- Fisher test including initially (Aug 2014) even colonies. Both colonies placed in the D-dominated category. P-value= 0.1064
Input =("
Year C D
2014 15 14
2016 8 21
")
Matriz <- as.matrix(read.table(textConnection(Input),
header=TRUE,
row.names=1))
fisher.test(Matriz,
alternative="two.sided")##
## Fisher's Exact Test for Count Data
##
## data: Matriz
## p-value = 0.1064
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.8336155 9.7658946
## sample estimates:
## odds ratio
## 2.760908- Fisher test including initially (Aug 2014) even colonies. One colony was placed in the C-dominated category (D proportion=0.47), while the second one was placed in the D-dominated category (D proportion = 0.50). P-value = 0.06099
Input =("
Year C D
2014 16 13
2016 8 21
")
Matriz <- as.matrix(read.table(textConnection(Input),
header=TRUE,
row.names=1))
fisher.test(Matriz,
alternative="two.sided")##
## Fisher's Exact Test for Count Data
##
## data: Matriz
## p-value = 0.06099
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.9556883 11.2396870
## sample estimates:
## odds ratio
## 3.162767- Fisher test removing initially (Aug 2014) even colonies. N=27 colonies instead 29. P-value = 0.09778
Input =("
Year C D
2014 15 12
2016 8 19
")
Matriz <- as.matrix(read.table(textConnection(Input),
header=TRUE,
row.names=1))
fisher.test(Matriz,
alternative="two.sided")##
## Fisher's Exact Test for Count Data
##
## data: Matriz
## p-value = 0.09778
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.848568 10.659493
## sample estimates:
## odds ratio
## 2.906792- Fisher test including initially (Aug 2014) even colonies. Both colonies placed in the C-dominated category (D proportion=0.47, colony 267; and 0.50, colony 269). p-value = 0.03295
Input =("
Year C D
2014 17 12
2016 8 21
")
Matriz <- as.matrix(read.table(textConnection(Input),
header=TRUE,
row.names=1))
fisher.test(Matriz,
alternative="two.sided")##
## Fisher's Exact Test for Count Data
##
## data: Matriz
## p-value = 0.03295
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 1.094473 13.006776
## sample estimates:
## odds ratio
## 3.6297776a. Changes in colony proportions based on detected genera (qPCR)
- Based on detected Symbiodiniacea genera in 2014 (qPCR)
Histo_Dprop <- ggplot(data.Uva, aes (D.Prp , fill=(D14Community))) + ggthe_bw +
geom_histogram(binwidth = 0.05, boundary=0) +
facet_grid(~Year.2, labeller = labeller(Year.2=c(
`2014` = "August 2014", `2015` = "August 2015", `2016` = "April 2016"))) +
scale_fill_manual(values=c("blue3","purple", "red3"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "Initial symbiont community (August 2014)",
labels=c("Cladocopium-only", "Mixed", "Durusdinium-only")) +
scale_x_continuous(name= expression(
italic(Durusdinium)~proportion~(D/H)/(S/H)),
breaks=(seq(0,1,by=0.1))) +
scale_y_continuous((name= "Number of colonies"),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
Histo_Dprop
Traject_Proportion_Facet <- ggplot (data.Uva) +
ggthe_bw +
geom_line(aes (x=Year.2, D.Prp, colour=factor(Colony)))+
facet_wrap(~D14Community, labeller = labeller(
D14Community=c(`Cladocopium-only` = "Cladocopium-only (2014)",
`Mixed` = "Mixed (2014)",
`Durusdinium-only` = "Durusdinium-only (2014)"))) +
geom_jitter(aes(Year.2, D.Prp, colour=factor(Colony)), alpha=0.3,
height = 0.05, width = 0.15) +
#labs(x = "Año del muestreo", y = "Durusdinium Proportion") +
theme(legend.position="none")
Traject_Proportion_Facet
Figure Not used: qPCR and colony proportion by dominant community
Histo2014_Dprop <- ggplot(data.Uva[which(data.Uva$Year.2==2014), ], aes (D.Prp , fill=(D14Community))) +
ggthe_bw + coord_flip()+ ggtitle("a.")+
geom_histogram(binwidth = 0.05, boundary=0) +
scale_fill_manual(values=c("blue3","purple", "red3"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "Initial symbiont community (August 2014)",
labels=c("Cladocopium-only", "Mixed", "Durusdinium-only")) +
scale_x_continuous(name= expression(
italic(Durusdinium)~proportion~((D/H)/(S/H))),
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies \n (Aug 2014)"),
limits = c(0,17),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
#Histo2014_Dprop
Histo2016_Dprop <- ggplot(data.Uva[which(data.Uva$Year.2==2016), ], aes (D.Prp , fill=(D14Community))) +
ggthe_bw + coord_flip()+ ggtitle("c.")+
geom_histogram(binwidth = 0.05, boundary=0) +
scale_fill_manual(values=c("blue3","purple", "red3"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "Initial symbiont community (August 2014)",
labels=c("Cladocopium-only", "Mixed", "Durusdinium-only")) +
scale_x_continuous(name="",
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies \n (Aug 2016)"),
limits = c(0,17),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
#Histo2016_Dprop
Traject_Proportion<- ggplot(data.Uva, aes (x=Year.2, y=D.Prp , colour=(D14Community), group=(Colony))) + ggthe_bw + ggtitle("b.") +
geom_line(linetype=2) +
geom_jitter(height = 0.01, width = 0.25, alpha=0.4) +
scale_colour_manual(values=c("blue3","purple", "red3"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "Initial symbiont community (August 2014)",
labels=c("Cladocopium-only", "Mixed", "Durusdinium-only")) +
scale_x_continuous(name="",
breaks = c(2014, 2015, 2016),
labels = c(" Aug \n 2014", " Aug \n 2015", " Apr \n 2016")) +
scale_y_continuous((name= ""),
breaks=seq(0,1,by=0.1),
expand=c(0.01,0.01)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
#Traject_Proportion
Figure_2<-grid.arrange(Histo2014_Dprop,Traject_Proportion, Histo2016_Dprop,
nrow=1, widths=c(3/9, 3/9, 3/9))
# ggsave(file="Outputs/Figure_2qPCR.svg", plot=Figure_2, width=7, height=4)7a. Changes in total Symbiodiniaceae abundance (bleaching): qPCR
Based on initial (2014) community
SH_2014Co<-lmerTest::lmer(logTot.SH ~ Year * D14Community + (1|Colony), data=data.Uva)
summary(SH_2014Co)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: logTot.SH ~ Year * D14Community + (1 | Colony)
## Data: data.Uva
##
## REML criterion at convergence: 83.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6538 -0.5202 -0.1347 0.3700 2.9430
##
## Random effects:
## Groups Name Variance Std.Dev.
## Colony (Intercept) 0.03036 0.1742
## Residual 0.10764 0.3281
## Number of obs: 87, groups: Colony, 29
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) -1.6402 0.1238 71.1153 -13.246
## Year15 -0.4132 0.1547 52.0000 -2.672
## Year16 -0.5536 0.1547 52.0000 -3.579
## D14CommunityMixed -0.0137 0.1638 71.1153 -0.084
## D14CommunityDurusdinium-only 0.1232 0.1805 71.1153 0.682
## Year15:D14CommunityMixed 0.7731 0.2046 52.0000 3.778
## Year16:D14CommunityMixed 0.6960 0.2046 52.0000 3.402
## Year15:D14CommunityDurusdinium-only 0.4584 0.2255 52.0000 2.033
## Year16:D14CommunityDurusdinium-only 0.3869 0.2255 52.0000 1.716
## Pr(>|t|)
## (Intercept) < 2e-16 ***
## Year15 0.010053 *
## Year16 0.000756 ***
## D14CommunityMixed 0.933599
## D14CommunityDurusdinium-only 0.497283
## Year15:D14CommunityMixed 0.000408 ***
## Year16:D14CommunityMixed 0.001294 **
## Year15:D14CommunityDurusdinium-only 0.047152 *
## Year16:D14CommunityDurusdinium-only 0.092125 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Year15 Year16 D14CmM D14CD- Y15:D14CM Y16:D14CM Y15:D14CD
## Year15 -0.624
## Year16 -0.624 0.500
## D14CmmntyMx -0.756 0.472 0.472
## D14CmmntyD- -0.686 0.428 0.428 0.519
## Yr15:D14CmM 0.472 -0.756 -0.378 -0.624 -0.324
## Yr16:D14CmM 0.472 -0.378 -0.756 -0.624 -0.324 0.500
## Yr15:D14CD- 0.428 -0.686 -0.343 -0.324 -0.624 0.519 0.259
## Yr16:D14CD- 0.428 -0.343 -0.686 -0.324 -0.624 0.259 0.519 0.500 #step(SH_2014Co)
SH.emmc<-emmeans(SH_2014Co, ~Year*D14Community)
SH_groups<-multcomp::cld(SH.emmc)
SH_groups<-SH_groups[order(
SH_groups$Year, SH_groups$D14Community),]
SH_groups plot(emmeans(SH_2014Co, ~Year|~D14Community), comparisons = TRUE) +
coord_flip(xlim = NULL, ylim = NULL, expand = TRUE) + facet_wrap(~D14Community)
SH_2014Com <- ggplot(data.Uva, aes (Year.2, logTot.SH, colour=factor(D14Community))) + ggthe_bw +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 1,
position=position_dodge(width=0.95))+
stat_summary(fun.data = "mean_cl_boot",geom = "point",
position=position_dodge(width=0.95))+
stat_summary(fun.y=mean, geom="line", position=position_dodge(width=0.95)) +
scale_colour_manual(values=c("blue","purple", "red"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "Symbionts detectected in 2014",
labels=c("Cladocopium", "Mixed", "Durusdinium"))+
theme(legend.position=c(0.5,0.1),
strip.background =element_rect(fill="white")) +
scale_y_continuous(limits = c(-3, -0.8),
name="\n Total symbiont to host cell ratio (S/H)",
labels = (math_format(10^.x)),
expand = c(0.0, 0)) +
scale_x_discrete(name="Year")
SH_2014Com
Based on current community
SH_Community<-lmerTest::lmer(logTot.SH ~ Year * Community + (1|Colony), data=data.Uva)
summary(SH_Community)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: logTot.SH ~ Year * Community + (1 | Colony)
## Data: data.Uva
##
## REML criterion at convergence: 91.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.3973 -0.5171 -0.0987 0.3905 4.0587
##
## Random effects:
## Groups Name Variance Std.Dev.
## Colony (Intercept) 0.02404 0.1550
## Residual 0.12526 0.3539
## Number of obs: 87, groups: Colony, 29
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -1.60009 0.12804 76.71004 -12.497 < 2e-16
## Year15 -0.20778 0.15540 49.73821 -1.337 0.187280
## Year16 -0.71036 0.17979 47.25498 -3.951 0.000258
## CommunityMixed -0.09453 0.16829 77.96039 -0.562 0.575910
## CommunityDurusdinium-only 0.09882 0.18645 77.17622 0.530 0.597641
## Year15:CommunityMixed 0.58170 0.24003 54.72489 2.423 0.018708
## Year16:CommunityMixed 0.81742 0.24650 54.96060 3.316 0.001622
## Year15:CommunityDurusdinium-only 0.27313 0.22958 47.93526 1.190 0.240035
## Year16:CommunityDurusdinium-only 0.60594 0.24059 49.84010 2.519 0.015039
##
## (Intercept) ***
## Year15
## Year16 ***
## CommunityMixed
## CommunityDurusdinium-only
## Year15:CommunityMixed *
## Year16:CommunityMixed **
## Year15:CommunityDurusdinium-only
## Year16:CommunityDurusdinium-only *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Year15 Year16 CmmntM CmmnD- Y15:CM Y16:CM Y15:CD
## Year15 -0.733
## Year16 -0.602 0.497
## CommuntyMxd -0.754 0.576 0.457
## CmmntyDrsd- -0.686 0.504 0.414 0.524
## Yr15:CmmntM 0.472 -0.658 -0.322 -0.628 -0.326
## Yr16:CmmntM 0.459 -0.365 -0.740 -0.608 -0.318 0.412
## Yr15:CmmnD- 0.496 -0.677 -0.337 -0.383 -0.717 0.440 0.243
## Yr16:CmmnD- 0.451 -0.369 -0.747 -0.326 -0.682 0.240 0.543 0.541 #step(SH_Community)
SH.emmc<-emmeans(SH_Community, ~Year*Community)
SH_groups<-multcomp::cld(SH.emmc)
SH_groups<-SH_groups[order(
SH_groups$Year, SH_groups$Community),]
SH_groups #write.csv(SH_groups,"SH_groupsB.csv")
plot(emmeans(SH_Community, ~Year|~Community), comparisons = TRUE) +
coord_flip(xlim = NULL, ylim = NULL, expand = TRUE) + facet_wrap(~Community)
SH_Com <- ggplot(data.Uva, aes (Year.2, logTot.SH, colour=factor(Community))) + ggthe_bw +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 1,
position=position_dodge(width=0.95))+
stat_summary(fun.data = "mean_cl_boot",geom = "point",
position=position_dodge(width=0.95))+
stat_summary(fun.y=mean, geom="line", position=position_dodge(width=0.95)) +
scale_colour_manual(values=c("blue","purple", "red"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "Symbionts detectec in 2014",
labels=c("Cladocopium", "Mixed", "Durusdinium"))+
theme(legend.position=c(0.5,0.1),
strip.background =element_rect(fill="white")) +
scale_y_continuous(limits = c(-3, -0.8),
name="\n Total symbiont to host cell ratio (S/H)",
labels = (math_format(10^.x)),
expand = c(0.0, 0)) +
scale_x_discrete(name="Year")
SH_Com
anova(SH_2014Co, SH_Community)Considering initial and current symbiont community: Figure 3: qPCR
SH_Communities<-lmerTest::lmer(logTot.SH ~ Year * D14Community * Community + (1|Colony), data=data.Uva)
summary(SH_Communities)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: logTot.SH ~ Year * D14Community * Community + (1 | Colony)
## Data: data.Uva
##
## REML criterion at convergence: 81.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.73407 -0.48463 -0.06579 0.35339 2.91935
##
## Random effects:
## Groups Name Variance Std.Dev.
## Colony (Intercept) 0.02179 0.1476
## Residual 0.11026 0.3320
## Number of obs: 87, groups: Colony, 29
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) -1.6402 0.1211 70.1550 -13.541
## Year15 -0.4132 0.1565 47.2363 -2.640
## Year16 -0.7093 0.1688 50.6884 -4.202
## D14CommunityMixed -0.7146 0.3265 73.0936 -2.189
## D14CommunityDurusdinium-only -0.6296 0.3928 72.9975 -1.603
## CommunityMixed 0.7009 0.2845 70.5547 2.464
## CommunityDurusdinium-only 0.7527 0.3508 70.9827 2.146
## Year15:D14CommunityMixed 1.5019 0.3806 64.8657 3.946
## Year16:D14CommunityMixed 0.8259 0.2396 53.9251 3.447
## Year15:D14CommunityDurusdinium-only 1.2674 0.5110 61.4046 2.480
## Year16:D14CommunityDurusdinium-only 0.5426 0.2368 48.9752 2.292
## Year15:CommunityMixed -0.7379 0.3653 70.6556 -2.020
## Year15:CommunityDurusdinium-only -0.8090 0.4572 64.7392 -1.770
## Pr(>|t|)
## (Intercept) < 2e-16 ***
## Year15 0.011206 *
## Year16 0.000107 ***
## D14CommunityMixed 0.031805 *
## D14CommunityDurusdinium-only 0.113266
## CommunityMixed 0.016178 *
## CommunityDurusdinium-only 0.035332 *
## Year15:D14CommunityMixed 0.000198 ***
## Year16:D14CommunityMixed 0.001108 **
## Year15:D14CommunityDurusdinium-only 0.015876 *
## Year16:D14CommunityDurusdinium-only 0.026258 *
## Year15:CommunityMixed 0.047162 *
## Year15:CommunityDurusdinium-only 0.081500 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## fit warnings:
## fixed-effect model matrix is rank deficient so dropping 14 columns / coefficients #step(SH_Communities)
SH.emmc<-emmeans(SH_Communities, ~Year*D14Community*Community)
SH_groups<-cld(SH.emmc)
SH_groups<-as.data.frame(SH_groups[complete.cases(SH_groups),])
SH_groups<-SH_groups[order(SH_groups$Year,
SH_groups$D14Community,
SH_groups$Community),]
SH_groups #write.csv(SH_groups,"SH_Communities_multicomp.csv")
SH_Comties <- ggplot(data.Uva, aes (Year.2, logTot.SH, colour=factor(Community))) + ggthe_bw+
facet_grid(~D14Community, labeller = labeller(D14Community=c(
`Cladocopium-only` = "Cladocopium-only (2014)",
`Mixed` = "Mixed (2014)",
`Durusdinium-only` = "Durusdinium-only (2014)"))) +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2, position=position_dodge(width=0.3))+
stat_summary(fun.data = "mean_cl_boot",geom = "point", position=position_dodge(width=0.3), alpha=0.5)+
stat_summary(fun.y=mean, geom="line", position=position_dodge(width=0.3)) +
scale_colour_manual(values=c("blue3","purple", "red3"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Detected symbionts",
labels=c("Cladocopium", "Mixed", "Durusdinium"))+
theme(legend.position=c(0.82,0.2),
legend.box.background = element_rect(),
strip.background =element_rect(fill="white"))+
scale_y_continuous(limits = c(-3, -0.8),
name="\n Total symbiont to host cell ratio (S/H)",
labels = (math_format(10^.x)),
expand = c(0.0, 0)) +
scale_x_continuous(name="",
breaks = c(14, 15, 16),
labels = c(" Aug \n 2014", " Aug \n 2015", " Apr \n 2016"))
SH_Comties
#ggsave(file="Figure_3.svg", plot=SH_Comties, width=6, height=4)8a. Changes in Cladocopium abundance
CH_Com <- ggplot(data.Uva, aes (Year.2, logC.SH, colour=factor(DominantCommunity))) + ggthe_bw+ ggtitle("a.")+
#geom_point()+
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2, position=position_dodge(width=0.3))+
stat_summary(fun.data = "mean_cl_boot",geom = "point", position=position_dodge(width=0.3), alpha=0.5)+
stat_summary(fun.y=mean, geom="line", position=position_dodge(width=0.3)) +
scale_colour_manual(values=c("blue3","purple", "red3"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Detected symbionts",
labels=c("Cladocopium", "Mixed", "Durusdinium"))+
theme(legend.position=c(0.82,0.8),
legend.box.background = element_rect(),
strip.background =element_rect(fill="white"))+
scale_y_continuous(limits = c(-3, -0.7),
name="\n Cladocopium to host cell ratio (C/H)",
labels = (math_format(10^.x)),
expand = c(0.0, 0)) +
scale_x_continuous(name="",
breaks = c(2014, 2015, 2016),
labels = c(" ", " ", " "))
#CH_Com
CH_Com2<-CH_Com + facet_grid(~D14Dominant, labeller = labeller(D14Dominant=c(
`Cladocopium-dominated` = "Cladocopium- \n dominated (2014)",
`Durusdinium-dominated` = "Durusdinium- \n dominated (2014)")))
CH<-lmer(logC.SH ~ Year* Community* D14Dominant + (1|Colony), data=data.Uva)
summary(CH)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: logC.SH ~ Year * Community * D14Dominant + (1 | Colony)
## Data: data.Uva
##
## REML criterion at convergence: 113.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.70311 -0.36865 -0.09556 0.50582 1.89250
##
## Random effects:
## Groups Name Variance Std.Dev.
## Colony (Intercept) 0.1182 0.3438
## Residual 0.4032 0.6350
## Number of obs: 55, groups: Colony, 21
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) -1.54068 0.23810 44.78323 -6.471
## Year15 -0.25625 0.28031 26.58185 -0.914
## Year16 -0.73614 0.32339 24.71137 -2.276
## CommunityMixed -0.32887 0.37022 44.99572 -0.888
## D14DominantDurusdinium-dominated -1.09559 0.41353 43.96389 -2.649
## Year15:CommunityMixed 0.59388 0.61658 32.70059 0.963
## Year16:CommunityMixed -0.03386 0.50241 28.81910 -0.067
## Year15:D14DominantDurusdinium-dominated -0.54917 0.68606 30.91184 -0.800
## Year16:D14DominantDurusdinium-dominated 0.88798 0.66184 31.43058 1.342
## Pr(>|t|)
## (Intercept) 6.36e-08 ***
## Year15 0.3688
## Year16 0.0318 *
## CommunityMixed 0.3791
## D14DominantDurusdinium-dominated 0.0112 *
## Year15:CommunityMixed 0.3425
## Year16:CommunityMixed 0.9467
## Year15:D14DominantDurusdinium-dominated 0.4296
## Year16:D14DominantDurusdinium-dominated 0.1893
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Year15 Year16 CmmntM D14DD- Y15:CM Y16:CM Y15:D1
## Year15 -0.715
## Year16 -0.579 0.492
## CommuntyMxd -0.622 0.502 0.369
## D14DmnntDr- -0.019 -0.038 0.003 -0.537
## Yr15:CmmntM 0.318 -0.474 -0.224 -0.512 0.275
## Yr16:CmmntM 0.400 -0.331 -0.659 -0.643 0.345 0.357
## Yr15:D14DD- 0.006 0.017 0.000 0.255 -0.468 -0.705 -0.185
## Yr16:D14DD- -0.021 0.011 0.011 0.308 -0.509 -0.162 -0.437 0.306
## fit warnings:
## fixed-effect model matrix is rank deficient so dropping 3 columns / coefficients #step(CH)
CH.emmc<-emmeans(CH, ~Year*Community* D14Dominant)
CH_groups<-cld(CH.emmc)
CH_groups<-as.data.frame(CH_groups[complete.cases(CH_groups),])
CH_groups<-CH_groups[order(CH_groups$Year,CH_groups$D14Dominant, CH_groups$Community),]
CH_groups #write.csv(CH_groups,"CH_groups.csv")9a. Changes in Durusdinium abundance
DH_Com <- ggplot(data.Uva, aes (Year.2, logD.SH, colour=factor(Community))) +
ggthe_bw+ ggtitle("b.")+
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2, position=position_dodge(width=0.3))+
stat_summary(fun.data = "mean_cl_boot",geom = "point", position=position_dodge(width=0.3), alpha=0.5)+
stat_summary(fun.y=mean, geom="line", position=position_dodge(width=0.3)) +
scale_colour_manual(values=c("purple", "red3"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Detected symbionts",
labels=c("Mixed", "Durusdinium"))+
theme(legend.position=c(0.82,0.2),
legend.box.background = element_rect(),
strip.background =element_rect(fill="white"))+
scale_y_continuous(limits = c(-6, -0.7),
name="\n Durusdinium to host cell ratio (D/H)",
labels = (math_format(10^.x)),
expand = c(0.0, 0)) +
scale_x_continuous(name="",
breaks = c(2014, 2015, 2016),
labels = c(" Aug \n 2014", " Aug \n 2015", " Apr \n 2016"))
# DH_Com
DH_Com2<-DH_Com + facet_grid(~D14Dominant, labeller = labeller(D14Dominant=c(
`Cladocopium-dominated` = "Cladocopium- \n dominated (2014)",
`Durusdinium-dominated` = "Durusdinium- \n dominated (2014)")))
DH_LM<-lmer(logD.SH ~ Year* Community* D14Dominant + (1|Colony), data=data.Uva)
summary(DH_LM)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: logD.SH ~ Year * Community * D14Dominant + (1 | Colony)
## Data: data.Uva
##
## REML criterion at convergence: 66.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.0953 -0.5756 0.0079 0.3844 3.6310
##
## Random effects:
## Groups Name Variance Std.Dev.
## Colony (Intercept) 0.0000 0.0000
## Residual 0.1687 0.4107
## Number of obs: 58, groups: Colony, 22
##
## Fixed effects:
## Estimate Std. Error
## (Intercept) -4.86573 0.16766
## Year15 0.07049 0.33532
## Year16 3.21459 0.23711
## CommunityDurusdinium-only 0.27806 0.50992
## D14DominantDurusdinium-dominated 3.15614 0.23711
## Year15:CommunityDurusdinium-only -0.35506 0.32897
## Year16:CommunityDurusdinium-only -0.22094 0.38416
## Year15:D14DominantDurusdinium-dominated 0.36688 0.42745
## Year16:D14DominantDurusdinium-dominated -3.09697 0.41068
## CommunityDurusdinium-only:D14DominantDurusdinium-dominated -0.08557 0.45916
## df t value
## (Intercept) 48.00000 -29.021
## Year15 48.00000 0.210
## Year16 48.00000 13.557
## CommunityDurusdinium-only 48.00000 0.545
## D14DominantDurusdinium-dominated 48.00000 13.311
## Year15:CommunityDurusdinium-only 48.00000 -1.079
## Year16:CommunityDurusdinium-only 48.00000 -0.575
## Year15:D14DominantDurusdinium-dominated 48.00000 0.858
## Year16:D14DominantDurusdinium-dominated 48.00000 -7.541
## CommunityDurusdinium-only:D14DominantDurusdinium-dominated 48.00000 -0.186
## Pr(>|t|)
## (Intercept) < 2e-16 ***
## Year15 0.834
## Year16 < 2e-16 ***
## CommunityDurusdinium-only 0.588
## D14DominantDurusdinium-dominated < 2e-16 ***
## Year15:CommunityDurusdinium-only 0.286
## Year16:CommunityDurusdinium-only 0.568
## Year15:D14DominantDurusdinium-dominated 0.395
## Year16:D14DominantDurusdinium-dominated 1.09e-09 ***
## CommunityDurusdinium-only:D14DominantDurusdinium-dominated 0.853
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Year15 Year16 CmmnD- D14DD- Y15:CD Y16:CD Y15:D1 Y16:D1
## Year15 -0.500
## Year16 -0.707 0.354
## CmmntyDrsd- 0.000 0.000 -0.232
## D14DmnntDr- -0.707 0.354 0.500 -0.232
## Yr15:CmmnD- 0.000 0.000 0.000 -0.293 0.360
## Yr16:CmmnD- 0.000 0.000 0.000 -0.753 0.309 0.389
## Yr15:D14DD- 0.392 -0.784 -0.277 0.129 -0.555 -0.500 -0.171
## Yr16:D14DD- 0.408 -0.204 -0.577 0.671 -0.577 -0.208 -0.713 0.320
## CmD-:D14DD- 0.000 0.000 0.258 -0.900 0.000 0.000 0.558 0.000 -0.596
## fit warnings:
## fixed-effect model matrix is rank deficient so dropping 2 columns / coefficients
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular #step(DH_LM)
DH.emmc<-emmeans(DH_LM, ~Year*Community* D14Dominant)
DH_groups<-cld(DH.emmc)
DH_groups<-as.data.frame(DH_groups[complete.cases(DH_groups),])
DH_groups<-DH_groups[order(DH_groups$Year,DH_groups$D14Dominant, DH_groups$Community),]
DH_groups #write.csv(DH_groups,"DH_groups.csv")Figure 4: Changes in C/H and D/H
Figure_4<-grid.arrange(CH_Com2, DH_Com2,nrow=2)
# ggsave(file="Outputs/Figure_4.svg", plot=Figure_4, width=7, height=4)Symbiodiiaceae ITS2 type + ORF lineage data analysis
Added Pocillopora ORF lineage and Symbiodiniaceae ITS2 type hosted in 2014
#data.Uva$Counts<-paste(data.Uva$DominantCommunity, data.Uva$ORFCommunity, sep = "_")
data.Uva$ORFCommunity<-paste(data.Uva$ORF, data.Uva$Community, sep = "_")
data.Uva$ORFCommunity<-factor(data.Uva$ORFCommunity,
levels = c("Type_3_Cladocopium-only",
"Type_1_Cladocopium-only",
"Type_1_Mixed",
"Type-1_Durusdinium-only"))
data.Uva$ORFCommunity2<-paste(data.Uva$ORF,
data.Uva$DominantCommunity, sep = "_")
data.Uva$ORFCommunity2<-as.factor(data.Uva$ORFCommunity2)
data.Uva$DominantITS2[which(
(data.Uva$DominantCommunity=="Durusdinium-dominated") &
(data.Uva$ORF=="Type_1")
)] <- "D1-dominated"
data.Uva$DominantITS2[which(
(data.Uva$DominantCommunity=="Cladocopium-dominated") &
(data.Uva$ORF=="Type_1")
)] <- "C1bc-dominated"
data.Uva$DominantITS2[which(
(data.Uva$DominantCommunity=="Cladocopium-dominated") &
(data.Uva$ORF=="Type_3")
)] <- "C1d-dominated"
data.Uva$DominantITS2<-as.factor(data.Uva$DominantITS2)
Summary_Community4<-data.Uva %>%
group_by(Year.2, DominantITS2) %>%
dplyr::summarise(total.count=n())
Summary_Community4 Summary_Community5<-data.Uva %>%
group_by(ORF, Year.2, DominantITS2) %>%
dplyr::summarise(total.count=n())
Summary_Community55b. Changes in the number of colonies dominated by each genus
- Fisher test including initially (Aug 2014) 50/50 colonies. Both colonies placed in the D-dominated category to be conservative.
# Only type 1 colonies (n=23), p-value = 0.03514
Input =("
Year C D
2014 9 14
2016 2 21
")
Matriz <- as.matrix(read.table(textConnection(Input),
header=TRUE,
row.names=1))
fisher.test(Matriz,
alternative="two.sided")##
## Fisher's Exact Test for Count Data
##
## data: Matriz
## p-value = 0.03514
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 1.105302 70.563784
## sample estimates:
## odds ratio
## 6.4764# Type 1 and type 3 colonies: p-value = 0.1064
Input =("
Year C D
2014 15 14
2016 8 21
")
Matriz <- as.matrix(read.table(textConnection(Input),
header=TRUE,
row.names=1))
fisher.test(Matriz,
alternative="two.sided")##
## Fisher's Exact Test for Count Data
##
## data: Matriz
## p-value = 0.1064
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.8336155 9.7658946
## sample estimates:
## odds ratio
## 2.7609086b. Changes in colony proportions based on ITS2 profiles (2014) and genus dominance
- Histograms based on ITS2 detected symbionts in 2014
data.Uva$ITS2_14<-factor(data.Uva$ITS2_14, levels = c("C1d", "C1b-c", "D1"))
data.Uva$ITS2_16<-factor(data.Uva$ITS2_16, levels = c("C1d", "C1b-c", "D1"))
Histo_Dprop_ITS2 <- ggplot(data.Uva, aes (D.Prp , fill=(ITS2_14)), shape=(D14Community)) + ggthe_bw +
geom_histogram(binwidth = 0.05, boundary=0) +
facet_grid(~Year.2, labeller = labeller(Year.2=c(
`2014` = "August 2014", `2015` = "August 2015", `2016` = "April 2016"))) +
scale_fill_manual(values=c("blue3", "purple", "red3"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "Initial ITS2 community (August 2014)",
labels=c("Type 3 + C1d", "Type 1 + C1b-c", "Type 1 + D1")) +
scale_x_continuous(name= expression(
italic(Durusdinium)~proportion~(D/H)/(S/H)),
breaks=(seq(0,1,by=0.1))) +
scale_y_continuous((name= "Number of colonies"),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
Histo_Dprop_ITS2 +coord_flip()
- Histograms based on ORF lineage
Histo_Dprop_ORF <- ggplot(data.Uva, aes (D.Prp , fill=(ORF)), shape=(D14Community)) + ggthe_bw +
geom_histogram(binwidth = 0.05, boundary=0) +
facet_grid(~Year.2, labeller = labeller(Year.2=c(
`2014` = "August 2014",
`2015` = "August 2015", `2016` = "April 2016"))) +
scale_fill_manual(values=c("black","gray"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "ORF type",
labels=c("Type 3", "Type 1")) +
scale_x_continuous(name= expression(
italic(Durusdinium)~proportion~(D/H)/(S/H)),
breaks=(seq(0,1,by=0.1))) +
scale_y_continuous((name= "Number of colonies"),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
Histo_Dprop_ORF +coord_flip()
- Durusdinium proportion trajectory
# By ITS2 type
Traject_Proportion_FacetITS2 <- ggplot (data.Uva) +
ggthe_bw +
geom_line(aes (x=Year.2, D.Prp, colour=factor(Colony)))+
facet_wrap(~ITS2_14) +
geom_jitter(aes(Year.2, D.Prp, colour=factor(Colony)), alpha=0.3,
height = 0.05, width = 0.15) +
#labs(x = "Año del muestreo", y = "Durusdinium Proportion") +
theme(legend.position="none")
Traject_Proportion_FacetITS2
# By Depth
Traject_Proportion_FacetORF <- ggplot (data.Uva) + ggthe_bw +
geom_line(aes (x=Year.2, D.Prp, group=Colony, linetype=(ORF)))+
geom_jitter(aes(Year.2, D.Prp, colour=(Depht), shape=ORF),
height = 0.00, width = 0.25, size=3, alpha=0.7) +
scale_colour_gradient(low = "red", high = "blue")
Traject_Proportion_FacetORF
Figure XX: ITS2
Histo2014_DpropITS2 <- ggplot(data.Uva[which(data.Uva$Year.2==2014), ], aes (D.Prp , fill=(ITS2_14))) +
ggthe_bw + ggtitle("a.")+
geom_histogram(binwidth = 0.05, boundary=0) +
scale_fill_manual(values=c("blue3","purple", "red3"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "Dominant ITS2 community (August 2014)",
labels=c("C1d", "C1b-c", "D.glynnii")) +
scale_x_continuous(name= expression(
italic(Durusdinium)~proportion),
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies \n (Aug 2014)"),
limits = c(0,17),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
#Histo2014_DpropITS2
Histo2016_DpropITS2 <- ggplot(data.Uva[which(data.Uva$Year.2==2016), ], aes (D.Prp , fill=(ITS2_16))) +
ggthe_bw + ggtitle("c.")+
geom_histogram(binwidth = 0.05, boundary=0) +
scale_fill_manual(values=c("blue3","purple", "red3"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "Dominant ITS2 community (April 2016)",
labels=c("C1d", "C1b-c", "D.glynnii")) +
scale_x_continuous(name="",
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies \n (Aug 2016)"),
limits = c(0,17),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
#Histo2016_DpropITS2
Traject_ProportionITS2<- ggplot(data.Uva, aes (x=Year.2, y=D.Prp ,
fill=ITS2_14, shape=Mortality, group=(Colony)))+
ggthe_bw + ggtitle("b.") +
geom_line(aes(colour=ITS2_14), linetype=2) +
geom_jitter(height = 0, width = 0.25, alpha=0.4) +
scale_colour_manual(values=c("blue3","purple", "red3"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "Initial ITS2 community (August 2014)",
labels=c("C1d", "C1b-c", "D.glynnii")) +
scale_fill_manual(values=c("blue3","purple", "red3"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "Initial ITS2 community (August 2014)",
labels=c("C1d", "C1b-c", "D.glynnii")) +
scale_shape_manual(values=c(21, 22, 24),
labels=c("No", "Partial", "Almost total"))+
scale_x_continuous(name="",
breaks = c(2014, 2015, 2016),
labels = c(" Aug \n 2014", " Aug \n 2015", " Apr \n 2016")) +
scale_y_continuous((name= ""),
breaks=seq(0,1,by=0.1),
expand=c(0.01,0.01)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
#Traject_ProportionITS2
Figure_2_ITS2<-grid.arrange(Histo2014_DpropITS2, Traject_ProportionITS2,Histo2016_DpropITS2,
nrow=1, widths=c(2/9, 5/9, 2/9))
# ggsave(file="Outputs/Figure_ITS2.svg", plot=Figure_2_ITS2, width=6.5, height=5)Histo2014_DpropITS2b<- Histo2014_DpropITS2 + facet_grid(ORF~.)
Histo2016_DpropITS2b<- Histo2016_DpropITS2 + facet_grid(ORF~.)
Traject_ProportionITS2b<- Traject_ProportionITS2 + facet_grid(ORF~.)
Figure_2_ITS2b<-grid.arrange(Histo2014_DpropITS2b, Traject_ProportionITS2b,Histo2016_DpropITS2b,
nrow=1, widths=c(3/9, 3/9, 3/9))
# ggsave(file="Outputs/Figure_2_ITS2b.svg", plot=Figure_2_ITS2b, width=7, height=4)Figure XX: ORF
Histo2014_DpropORF <- ggplot(data.Uva[which(data.Uva$Year.2==2014), ], aes (D.Prp , fill=(ORFCommunity))) +
ggthe_bw + coord_flip()+ ggtitle("a.")+
geom_histogram(binwidth = 0.05, boundary=0, colour="black") +
scale_fill_manual(values=c("gray", "white","snow", "mintcream"),
guide = guide_legend(title.position = "top", nrow = 2),
name = "")+
#labels=c("Type 3", "Type 1")) +
scale_x_continuous(#name= expression(
#italic(Durusdinium)~proportion~((D/H)/(S/H))),
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies \n (Aug 2014)"),
limits = c(0,16),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA),
panel.border = element_blank(),
axis.line = element_line(colour = "black"))
#Histo2014_DpropORF
Histo2016_DpropORF <- ggplot(data.Uva[which(data.Uva$Year.2==2016), ], aes (D.Prp , fill=(ORFCommunity))) +
ggthe_bw + coord_flip()+ ggtitle("c.")+
geom_histogram(binwidth = 0.05, boundary=0,
colour="black") +
scale_fill_manual(values=c("gray", "white","snow", "mintcream"),
guide = guide_legend(title.position = "top", nrow = 2),
name = "")+
#labels=c("Type 3", "Type 1")) +
scale_x_continuous(name="",
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies \n (Aug 2016)"),
limits = c(0,16),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA),
panel.border = element_blank(),
axis.line = element_line(colour = "black"))
#Histo2016_DpropORF
Traject_ProportionORF_a<- ggplot(data.Uva, aes (x=Year.2, y=D.Prp))+
ggthe_bw + ggtitle("b.") +
geom_line(aes(colour=(ORF), group=(Colony)), linetype=2) +
scale_colour_manual(values=c("black", "gray"),
guide = guide_legend(title.position = "top", nrow = 2),
name = "Pocillopora type",
labels=c("Type 3", "Type 1")) +
scale_x_continuous(name="",
breaks = c(2014, 2015, 2016),
labels = c(" Aug \n 2014", " Aug \n 2015", " Apr \n 2016")) +
scale_y_continuous((name= ""),
breaks=seq(0,1,by=0.1),
expand=c(0.01,0.01)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
Traject_ProportionORF <- Traject_ProportionORF_a +
geom_jitter((aes(colour=factor(ORF))),
height = 0, width = 0.3, color="black", size=1, shape=21) +
scale_fill_manual(values=c("white", "gray"))
#Traject_ProportionORF
Figure_2_ORF<-grid.arrange(Histo2014_DpropORF, Traject_ProportionORF, Histo2016_DpropORF,
nrow=1, widths=c(3/9, 3/9, 3/9))
# ggsave(file="Outputs/Figure_2_ORF.svg",plot=Figure_2_ORF, width=7, height=3.8, dpi=300)Figure XX: ITS2 + Moratality
data.Uva$Mor_ITS2<-paste(data.Uva$DominantITS2, data.Uva$Mortality, sep = "_")
data.Uva$Mor_ITS2<-factor(data.Uva$Mor_ITS2, levels = c(
"C1d-dominated_No", "C1d-dominated_Partial", "C1d-dominated_Almost total",
"C1bc-dominated_No", "C1bc-dominated_Partial",
"D1-dominated_No", "D1-dominated_Partial", "D1-dominated_Almost total"))
summary(data.Uva$Mor_ITS2)## C1d-dominated_No C1d-dominated_Partial
## 11 3
## C1d-dominated_Almost total C1bc-dominated_No
## 4 17
## C1bc-dominated_Partial D1-dominated_No
## 3 43
## D1-dominated_Partial D1-dominated_Almost total
## 3 3Histo2014_DpropITS2 <- ggplot(data.Uva[which(data.Uva$Year.2==2014), ], aes (D.Prp , fill=(Mor_ITS2))) +
ggthe_bw + coord_flip()+ ggtitle("a.")+
geom_histogram(binwidth = 0.05, boundary=0) +
scale_fill_manual(values=c("skyblue1",
"darkorchid1",
"red1"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Dominant ITS2 community (August 2014)") +
scale_x_continuous(name= expression(
italic(Durusdinium)~proportion~((D/H)/(S/H))),
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies \n (Aug 2014)"),
limits = c(0,17),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
#Histo2014_DpropITS2
Histo2016_DpropITS2 <- ggplot(data.Uva[which(data.Uva$Year.2==2016), ], aes (D.Prp , fill=(Mor_ITS2))) +
ggthe_bw + coord_flip()+ ggtitle("c.")+
geom_histogram(binwidth = 0.05, boundary=0) +
scale_fill_manual(values=c("skyblue1", "blue2", "blue4",
"darkorchid3",
"red1", "red3", "red4"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Dominant ITS2 community (April 2016)") +
scale_x_continuous(name="",
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies \n (Aug 2016)"),
limits = c(0,17),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
#Histo2016_DpropITS2
HistoAll_DpropITS2 <- ggplot(data.Uva, aes (D.Prp , fill=(Mor_ITS2))) +
ggthe_bw + coord_flip()+ ggtitle("c.")+
geom_histogram(binwidth = 0.05, boundary=0) +
scale_fill_manual(values=c("skyblue1", "blue2", "blue4",
"darkorchid1","darkorchid3",
"red1", "red3", "red4"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Dominant ITS2 community (April 2016)") +
scale_x_continuous(name="",
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies \n (Aug 2016)"),
limits = c(0,17),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA)) + facet_grid(ORF~Year)
#HistoAll_DpropITS2
Traject_ProportionITS2<- ggplot(data.Uva, aes (x=Year.2, y=D.Prp, group=(Colony)))+
ggthe_bw + ggtitle("b.") +
geom_line(aes(colour=ITS2_14),linetype=2) +
geom_jitter(aes(colour=Mor_ITS2, shape=Mortality),
height = 0.00, width = 0.25, alpha=0.9) +
scale_colour_manual(values=c(
"darkorchid1", "darkorchid1","darkorchid3",
"skyblue1","blue4", "skyblue1", "blue2",
"red1", "red4", "red1", "red3" ),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Dominant ITS2 community (April 2016)") +
scale_x_continuous(name="",
breaks = c(2014, 2015, 2016),
labels = c(" Aug \n 2014", " Aug \n 2015", " Apr \n 2016")) +
scale_y_continuous((name= ""),
breaks=seq(0,1,by=0.1),
expand=c(0.01,0.01)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))
#Traject_ProportionITS2
Figure_2_ITS2<-grid.arrange(Histo2014_DpropITS2, Traject_ProportionITS2,Histo2016_DpropITS2,
nrow=1, widths=c(3/9, 3/9, 3/9))
# ggsave(file="Outputs/Figure_ITS2_Mortality.svg", plot=Figure_2_ITS2, width=7, height=4)Figure 3b and 3c: ITS2 + Moratality
data.Uva$Mor_ITS2<-paste(data.Uva$DominantITS2, data.Uva$Mortality, sep = "_")
data.Uva$Mor_ITS2<-factor(data.Uva$Mor_ITS2, levels = c(
"C1d-dominated_Almost total", "C1d-dominated_Partial","C1d-dominated_No",
"C1bc-dominated_Partial","C1bc-dominated_No",
"D1-dominated_Almost total", "D1-dominated_Partial", "D1-dominated_No" ))
summary(data.Uva$Mor_ITS2)## C1d-dominated_Almost total C1d-dominated_Partial
## 4 3
## C1d-dominated_No C1bc-dominated_Partial
## 11 3
## C1bc-dominated_No D1-dominated_Almost total
## 17 3
## D1-dominated_Partial D1-dominated_No
## 3 43Histo2014_DpropITS2 <- ggplot(data.Uva[which(data.Uva$Year.2==2014), ], aes (D.Prp , fill=(Mor_ITS2))) +
ggthe_bw + ggtitle("a.")+
geom_histogram(binwidth = 0.05, boundary=0) +
scale_fill_manual(values=c("black",
"black",
"black"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Dominant ITS2 community (August 2014)") +
scale_x_continuous(name= expression(
italic(Durusdinium)~proportion~((D/H)/(S/H))),
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies \n (Aug 2014)"),
limits = c(0,17),
breaks=seq(0,15,by=5)) +
theme(legend.position="none",
legend.title.align=0.5,
strip.background =element_rect(fill=NA)) + facet_grid(ORF~.)
#Histo2014_DpropITS2
Histo2016_DpropITS2 <- ggplot(data.Uva[which(data.Uva$Year.2==2016), ], aes (D.Prp , fill=(Mor_ITS2))) +
ggthe_bw + ggtitle("c.")+
geom_histogram(binwidth = 0.05, boundary=0) +
scale_fill_manual(values=c("white", "gray", "black",
"gray",
"white", "gray", "black"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Dominant ITS2 community (April 2016)") +
scale_x_continuous(name="",
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies \n (Aug 2016)"),
limits = c(0,17),
breaks=seq(0,15,by=5)) +
theme(legend.position="none",
legend.title.align=0.5,
strip.background =element_rect(fill=NA)) + facet_grid(ORF~.)
#Histo2016_DpropITS2
HistoAll_DpropITS2 <- ggplot(data.Uva, aes (D.Prp , fill=(Mor_ITS2))) +
ggthe_bw + ggtitle("c.")+
geom_histogram(binwidth = 0.05, boundary=0) +
scale_fill_manual(values=c("skyblue1", "blue2", "blue4",
"darkorchid1","darkorchid3",
"red1", "red3", "red4"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Dominant ITS2 community (April 2016)") +
scale_x_continuous(name="",
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies \n (Aug 2016)"),
limits = c(0,17),
breaks=seq(0,15,by=5)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA)) + facet_grid(ORF~Year)
#HistoAll_DpropITS2
HistoAll_DpropITS2 <- ggplot(data.Uva, aes (D.Prp , fill=(Mor_ITS2))) +
ggthe_bw + ggtitle("c.")+
geom_histogram(binwidth = 0.05, boundary=0, color = "black") +
scale_fill_manual(values=c("white", "gray", "black",
"gray","black",
"white", "gray", "black"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Dominant ITS2 community (April 2016)") +
scale_x_continuous(name="",
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies"),
limits = c(0,14),
breaks=seq(0,14,by=2)) +
theme(legend.position="bottom",
legend.title.align=0.5,
strip.background =element_rect(fill=NA)) +
facet_grid(ORF~Year) + theme(legend.position = "none")
#HistoAll_DpropITS2
#ggsave(file="Outputs/Panels.svg", plot=HistoAll_DpropITS2, width=6, height=2.8)
Traject_ProportionITS2<- ggplot(data.Uva, aes (x=Year.2, y=D.Prp, group=(Colony)))+
ggthe_bw + ggtitle("b.") +
geom_line(linetype=2) +
geom_jitter(aes(colour=Mortality, shape=ORF), size=1.5,
height = 0.00, width = 0.25, alpha=1) +
scale_colour_manual(values=c("black", "gray","white")) +
scale_x_continuous(name="",
breaks = c(2014, 2015, 2016),
labels = c(" Aug \n 2014", " Aug \n 2015", " Apr \n 2016")) +
scale_y_continuous((name= ""),
breaks=seq(0,1,by=0.1),
expand=c(0.02,0.02)) +
theme(legend.position="none",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))# + facet_grid(ORF~.)
#Traject_ProportionITS2
#ggsave(file="Outputs/Trajectory.svg", plot=Traject_ProportionITS2, width=6.5, height=3)
Figure_2_ITS2<-grid.arrange(Histo2014_DpropITS2, Traject_ProportionITS2,Histo2016_DpropITS2,
nrow=1, widths=c(2.5/9, 4/9, 2.5/9))
# ggsave(file="Outputs/Figure_ITS2_Mortality2.svg", plot=Figure_2_ITS2, width=7, height=4)DpropORF <- ggplot(data.Uva, aes (D.Prp , fill=(Mortality))) +
ggthe_bw + ggtitle("a.")+
geom_histogram(binwidth = 0.1, boundary=0, color="black") +
scale_fill_manual(values=c("black",
"gray",
"white"),
guide = guide_legend(title.position = "top", nrow = 1),
name = "Colony condition") +
scale_x_continuous(name= expression(
italic(Durusdinium)~proportion~((D/H)/(S/H))),
breaks=(seq(0,1,by=0.1)),
expand=c(0.01,0.01)) +
scale_y_continuous((name= "Number of colonies
\n (Aug 2014)
\n (Aug 2014) \n (Aug 2014)"),
limits = c(0,17),
breaks=seq(0,15,by=5)) +
theme(legend.position="top",
legend.title.align=0.5,
strip.background =element_rect(fill=NA)) + facet_grid(ORF~Year)
#DpropORF
Traject_ORF<- ggplot(data.Uva, aes (x=Year.2, y=D.Prp, group=(Colony)))+
ggthe_bw + ggtitle("b.") +
geom_line(linetype=2) +
geom_jitter(aes(colour=Mortality, shape=ORF), size=2,
height = 0.00, width = 0.25, alpha=1) +
scale_colour_manual(values=c("black", "gray","white")) +
scale_x_continuous(name="",
breaks = c(2014, 2015, 2016),
labels = c(" Aug \n 2014", " Aug \n 2015", " Apr \n 2016")) +
scale_y_continuous((name= ""),
breaks=seq(0,1,by=0.1),
expand=c(0.02,0.02)) +
theme(legend.position="none",
legend.title.align=0.5,
strip.background =element_rect(fill=NA))# + facet_grid(ORF~.)
#Traject_ORF
Figure_3bc_ORF<-grid.arrange(Traject_ORF, DpropORF,
nrow=2, heights=c(1.2/3, 1.8/3))
#ggsave(file="Outputs/Fig3_v8.svg", plot=Figure_3bc_ORF, width=6.5, height=7)# Summary dominant genera
Community_prescence_absence<-data.Uva %>%
group_by(ORF, Community, Year.2) %>%
dplyr::summarise(total.count=n())
Community_prescence_absence Community_dominance<-data.Uva %>%
group_by(ORF, DominantCommunity, Year.2) %>%
dplyr::summarise(total.count=n())
Community_dominance Community_dominance<-data.Uva %>%
group_by(ORF, Year.2, DominantCommunity,
Community) %>%
dplyr::summarise(total.count=n())
Community_dominance # SH_mean<-data.Uva %>%
# group_by(ORF, Year.2, DominantCommunity) %>%
# dplyr::summarise(mean=mean(tot.SH))
# SH_mean
#
# SH_SE<-data.Uva %>%
# group_by(ORF, Year.2, DominantCommunity) %>%
# dplyr::summarise(sd=sd(tot.SH))
# SH_SE7b. Changes in total Symbiodiniaceae abundance (bleaching): qPCR + ORF
Figure 3: ORF + Dominant symbiont
#SH_2014CoI<-lm(logTot.SH ~ Year * ORF * DominantCommunity, data=data.Uva)
SH_2014CoI<-lmerTest::lmer(logTot.SH ~ Year * ORF * DominantCommunity + (1|Colony), data=data.Uva)
summary(SH_2014CoI)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: logTot.SH ~ Year * ORF * DominantCommunity + (1 | Colony)
## Data: data.Uva
##
## REML criterion at convergence: 75.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.3367 -0.5185 -0.0403 0.3713 3.4413
##
## Random effects:
## Groups Name Variance Std.Dev.
## Colony (Intercept) 0.01992 0.1411
## Residual 0.10475 0.3236
## Number of obs: 87, groups: Colony, 29
##
## Fixed effects:
## Estimate Std. Error df
## (Intercept) -1.68521 0.11691 77.39879
## Year15 0.22667 0.15257 51.72506
## Year16 -0.05117 0.26368 68.09623
## ORFType_3 0.02507 0.18560 75.84196
## DominantCommunityDurusdinium-dominated 0.14102 0.14919 77.99999
## Year15:ORFType_3 -0.79408 0.24123 51.72506
## Year16:ORFType_3 -0.72686 0.32318 62.79186
## Year15:DominantCommunityDurusdinium-dominated -0.06059 0.19555 51.72506
## Year16:DominantCommunityDurusdinium-dominated 0.01196 0.28860 68.27239
## t value Pr(>|t|)
## (Intercept) -14.415 <2e-16 ***
## Year15 1.486 0.1434
## Year16 -0.194 0.8467
## ORFType_3 0.135 0.8929
## DominantCommunityDurusdinium-dominated 0.945 0.3475
## Year15:ORFType_3 -3.292 0.0018 **
## Year16:ORFType_3 -2.249 0.0280 *
## Year15:DominantCommunityDurusdinium-dominated -0.310 0.7579
## Year16:DominantCommunityDurusdinium-dominated 0.041 0.9671
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Year15 Year16 ORFT_3 DmnCD- Y15:OR Y16:OR Y15:DC
## Year15 -0.653
## Year16 -0.373 0.289
## ORFType_3 -0.630 0.411 0.235
## DmnntCmmnD- -0.777 0.511 0.290 0.489
## Yr15:ORFT_3 0.413 -0.632 -0.183 -0.650 -0.323
## Yr16:ORFT_3 0.305 -0.236 -0.816 -0.483 -0.237 0.373
## Yr15:DmnCD- 0.509 -0.780 -0.226 -0.321 -0.655 0.493 0.184
## Yr16:DmnCD- 0.359 -0.264 -0.920 -0.226 -0.463 0.167 0.750 0.339
## fit warnings:
## fixed-effect model matrix is rank deficient so dropping 3 columns / coefficients anova(SH_2014CoI) ranova(SH_2014CoI) step(SH_2014CoI, reduce.random=FALSE)## Backward reduced random-effect table:
##
## Eliminated npar logLik AIC LRT Df Pr(>Chisq)
## <none> 11 -37.966 97.933
## (1 | Colony) 0 10 -38.884 97.768 1.8355 1 0.1755
##
## Backward reduced fixed-effect table:
## Degrees of freedom method: Satterthwaite
##
## Eliminated Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Year:ORF:DominantCommunity 0
##
## Model found:
## logTot.SH ~ Year * ORF * DominantCommunity + (1 | Colony) SH_2014CoI<-lm(logTot.SH ~ Year * ORF * DominantCommunity, data=data.Uva)
SH_I.emmc<-emmeans(SH_2014CoI, ~Year*ORF*DominantCommunity)
SH_I_groups<-multcomp::cld(SH_I.emmc)
SH_I_groups<-as.data.frame(SH_I_groups[complete.cases(SH_I_groups),])
SH_I_groups<-SH_I_groups[order(
SH_I_groups$Year, SH_I_groups$ORF),]
SH_I_groups #write.csv(SH_I_groups,"Outputs/SH_ORF_Dominant_groups2.csv")
# plot(emmeans(SH_2014CoI, ~Year|~DominantCommunity), comparisons = TRUE) + coord_flip(xlim = NULL, ylim = NULL, expand = TRUE) + facet_wrap(~DominantITS2)+ ggthe_bwSH_2014ComI <- ggplot(data.Uva, aes (Year.2, logTot.SH, colour=factor(DominantITS2))) + ggthe_bw + facet_wrap(~ORF) + ggtitle("a.")+
#geom_point(alpha=0.5)+
#geom_line(aes(group=(Colony)))+
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2,
position=position_dodge(width=0.5))+
stat_summary(fun.data = "mean_cl_boot",geom = "point",
position=position_dodge(width=0.5))+
stat_summary(fun.y=mean, geom="line", position=position_dodge(width=0.5)) +
scale_colour_manual(values=c("blue3","blue3" , "red3"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Dominant symbiont",
labels=c("Cladocopium", "Durusdinium"))+
theme(legend.position="mone",
strip.background =element_rect(fill="white"))+
scale_y_continuous(limits = c(-3, -0.1),
name="\n Total symbiont to host cell ratio (S/H)",
labels = (math_format(10^.x)),
expand = c(0.0, 0)) +
scale_x_continuous(name="",
breaks = c(2014, 2015, 2016),
labels = c(" ", "", " "))
SH_2014ComI
#ggsave(file="Figure_3.svg", plot=SH_ITS2_Comunities, width=6, height=4, dpi=300)8b. Changes in Cladocopium abundance
CH_Com2ITS2 <- ggplot(data.Uva, aes (Year.2, logC.SH, colour=factor(D14Dominant))) + ggthe_bw+ ggtitle("c.") + facet_grid(~ORF) + #geom_point(alpha=0.5)+
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2, position=position_dodge(width=0.3))+
stat_summary(fun.data = "mean_cl_boot",geom = "point", position=position_dodge(width=0.3))+
stat_summary(fun.y=mean, geom="line", position=position_dodge(width=0.3)) +
scale_colour_manual(values=c("blue3", "red3"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Dominant symbiont",
labels=c("Cladocopium", "Durusdinium"))+
theme(legend.position="mone",
strip.background =element_rect(fill="white"))+
scale_y_continuous(limits = c(-6, -0.1),
name="\n Cladocopium to host cell ratio (C/H)",
labels = (math_format(10^.x)),
expand = c(0.0, 0)) +
scale_x_continuous(name="",
breaks = c(2014, 2015, 2016),
labels = c(" ", " ", " "))
CH_Com2ITS2
CH_I<-lmerTest::lmer(logC.SH ~ Year* DominantITS2 + (1|Colony), data=data.Uva)
summary(CH_I)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: logC.SH ~ Year * DominantITS2 + (1 | Colony)
## Data: data.Uva
##
## REML criterion at convergence: 108.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.89915 -0.40700 -0.07625 0.38592 1.92091
##
## Random effects:
## Groups Name Variance Std.Dev.
## Colony (Intercept) 0.06891 0.2625
## Residual 0.38254 0.6185
## Number of obs: 55, groups: Colony, 21
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -1.68510 0.22354 44.93743 -7.538 1.65e-09
## Year15 0.22932 0.29156 28.78344 0.787 0.437988
## Year16 -0.12883 0.50336 39.26788 -0.256 0.799334
## DominantITS2C1d-dominated 0.02496 0.35386 44.57426 0.071 0.944077
## DominantITS2D1-dominated -1.27282 0.35244 45.62565 -3.611 0.000755
## Year15:DominantITS2C1d-dominated -0.79673 0.46100 28.78344 -1.728 0.094658
## Year16:DominantITS2C1d-dominated -0.64921 0.61716 35.78277 -1.052 0.299886
## Year15:DominantITS2D1-dominated -0.44588 0.49919 31.70568 -0.893 0.378478
## Year16:DominantITS2D1-dominated 0.31093 0.62644 42.90497 0.496 0.622184
##
## (Intercept) ***
## Year15
## Year16
## DominantITS2C1d-dominated
## DominantITS2D1-dominated ***
## Year15:DominantITS2C1d-dominated .
## Year16:DominantITS2C1d-dominated
## Year15:DominantITS2D1-dominated
## Year16:DominantITS2D1-dominated
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Year15 Year16 DITS2C DITS2D Y15:DITS2C Y16:DITS2C Y15:DITS2D
## Year15 -0.652
## Year16 -0.376 0.290
## DmnnITS2C1- -0.632 0.412 0.238
## DmnnITS2D1- -0.631 0.414 0.238 0.398
## Y15:DITS2C1 0.412 -0.632 -0.183 -0.651 -0.262
## Y16:DITS2C1 0.307 -0.236 -0.816 -0.486 -0.194 0.373
## Y15:DITS2D1 0.382 -0.584 -0.170 -0.241 -0.606 0.369 0.138
## Y16:DITS2D1 0.339 -0.233 -0.816 -0.214 -0.538 0.147 0.665 0.345 anova (CH_I) #step(CH_I)
CH_I.emmc<-emmeans(CH_I, ~Year*DominantITS2)
CH_I_groups<-cld(CH_I.emmc)
CH_I_groups<-as.data.frame(CH_I_groups[complete.cases(CH_I_groups),])
CH_I_groups<-CH_I_groups[order(CH_I_groups$DominantITS2,CH_I_groups$Year),]
CH_I_groups #write.csv(CH_I_groups,"CH_I_groups.csv")9b. Changes in Durusdinium abundance
DH_Com2ITS2 <- ggplot(data.Uva, aes (Year.2, logD.SH, colour=factor(D14Dominant))) + # geom_point(alpha=0.5)+
ggthe_bw+ ggtitle("d.")+ facet_grid(~ORF) +
stat_summary(fun.data = "mean_cl_boot",geom = "errorbar", width = 0.2, position=position_dodge(width=0.3))+
stat_summary(fun.data = "mean_cl_boot",geom = "point", position=position_dodge(width=0.3))+
stat_summary(fun.y=mean, geom="line", position=position_dodge(width=0.3)) +
scale_colour_manual(values=c("blue", "red3"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Dominant symbionts",
labels=c("Cladocopium-dominated", "Durusdinium-dominated"))+
theme(legend.position=c(0.75,0.5),
legend.box.background = element_rect(),
strip.background =element_rect(fill="white"))+
scale_y_continuous(limits = c(-6, -0.1),
name="\n Durusdinium to host cell ratio (D/H)",
labels = (math_format(10^.x)),
expand = c(0.0, 0)) +
scale_x_continuous(name="",
breaks = c(2014, 2015, 2016),
# labels = c(" Aug \n 2014", " Aug \n 2015", " Apr \n 2016"))
labels = c(" ", " ", " "))
DH_Com2ITS2
DH_LM_ITS2<-lmer(logD.SH ~ Year* DominantITS2 + (1|Colony), data=data.Uva)
summary(DH_LM_ITS2)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: logD.SH ~ Year * DominantITS2 + (1 | Colony)
## Data: data.Uva
##
## REML criterion at convergence: 62.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.1646 -0.5042 0.0345 0.3558 3.7511
##
## Random effects:
## Groups Name Variance Std.Dev.
## Colony (Intercept) 0.000 0.0000
## Residual 0.158 0.3975
## Number of obs: 58, groups: Colony, 22
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -4.86573 0.16229 52.00000 -29.981 < 2e-16
## Year15 0.07049 0.32459 52.00000 0.217 0.829
## Year16 2.95052 0.42939 52.00000 6.871 7.93e-09
## DominantITS2D1-dominated 3.26614 0.19398 52.00000 16.838 < 2e-16
## Year15:DominantITS2D1-dominated 0.14077 0.35768 52.00000 0.394 0.696
## Year16:DominantITS2D1-dominated -2.96086 0.45076 52.00000 -6.569 2.41e-08
##
## (Intercept) ***
## Year15
## Year16 ***
## DominantITS2D1-dominated ***
## Year15:DominantITS2D1-dominated
## Year16:DominantITS2D1-dominated ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Year15 Year16 DITS2D Y15:DI
## Year15 -0.500
## Year16 -0.378 0.189
## DmnnITS2D1- -0.837 0.418 0.316
## Y15:DITS2D1 0.454 -0.907 -0.171 -0.542
## Y16:DITS2D1 0.360 -0.180 -0.953 -0.430 0.233
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular #step(DH_LM_ITS2)
DH_I.emmc<-emmeans(DH_LM_ITS2, ~Year*DominantITS2)
DH_groups_I<-cld(DH_I.emmc)
DH_groups_I<-as.data.frame(DH_groups_I[complete.cases(DH_groups_I),])
DH_groups_I<-DH_groups_I[order(DH_groups_I$DominantITS2,DH_groups_I$Year),]
DH_groups_I #write.csv(DH_groups_I,"Output/DH_I_groups.csv")Figure 4: qPCR by Pocillopora lineage
Figure_4<-grid.arrange(SH_2014ComI, CH_Com2ITS2, DH_Com2ITS2,nrow=3)
Figure_4_InitailCommunity<-grid.arrange(CH_Com2ITS2, DH_Com2ITS2,nrow=2)
# ggsave(file="Outputs/Figure_3c.svg", plot=Figure_3, width=7, height=11)
# ggsave(file="Outputs/Figure_3d.png", plot=Figure_3_InitailCommunity, width=7, height=7)10. Supplmenmtary C/H D/H dynamics
Select data
log_Log.data<-data.Uva[,c("Colony","Year.2", "C.SH", "D.SH",
"D14Community", "Community","ITS2_14",
"ORF", "DominantITS2", "Mortality", "Mor_ITS2")]
log_Log.data<-log_Log.data %>%
filter(Year.2 != "2015")
log_Log.data$Yr<-as.factor(log_Log.data$Year.2)
log_Log.data$Yr_M<-as.factor(
paste(log_Log.data$Yr,
log_Log.data$Mortality,
sep = "_"))
log_Log.data$logC<-log10(log_Log.data$C.SH)
log_Log.data$logD<-log10(log_Log.data$D.SH)
log_Log.data$logC[which(log_Log.data$C.SH==0)] <- -6
log_Log.data$logD[which(log_Log.data$D.SH==0)] <- -6- ITS2-ORF
Community_Changes1 <- ggplot (log_Log.data, aes(logC, logD, gropup=Colony, colour=ITS2_14, shape=ORF)) +
geom_abline(intercept=0, slope=1, linetype=1, color="black")+
geom_abline(intercept=-1, slope=1, linetype=2, color="black")+
geom_abline(intercept=-2, slope=1, linetype=2, color="darkgray")+
geom_abline(intercept=-3, slope=1, linetype=2, color="gray")+
geom_abline(intercept=-4, slope=1, linetype=2, color="gray")+
geom_abline(intercept=-5, slope=1, linetype=2, color="gray")+
geom_abline(intercept=1, slope=1, linetype=2, color="black")+
geom_abline(intercept=2, slope=1, linetype=2, color="darkgray")+
geom_abline(intercept=3, slope=1, linetype=2, color="gray")+
geom_abline(intercept=4, slope=1, linetype=2, color="gray")+
geom_abline(intercept=5, slope=1, linetype=2, color="gray")+
scale_colour_manual(values=c("blue","purple", "red"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Dominant ITS2 types (2014)",
labels=c("C1d", "C1b-c", "D1"))+
ggthe_bw + theme(legend.position="bottom",
legend.box = "horizontal")+
scale_x_continuous(limits = c(-6,-0.1), "log10 (C/H cell ratio)")+
scale_y_continuous(limits = c(-6,-0.1), "log10 (D/H cell ratio)") +
geom_point(alpha=0.5, size=4) +
geom_path(size = 1,arrow =
arrow(type = "open", angle = 30, length = unit(0.1, "inches")))+
annotate("text", x = c(seq(-5.1, -1.1, by=1)), y = -0.9, label = c(
"10,0000D : 1C", "1,000D : 1C", "100D : 1C","10D : 1C", "C : D"), angle = 45)+
annotate("text", x = -0.9, y = c(seq(-4.7, -1.7, by=1)), label = c(
"1D : 10,000C", "1D : 1,000C", "1D : 100C", "1D : 10C"), angle = 45)
Community_Changes1
#ggsave(file="FigAndrewsFavoriteGraph_ITS2.png", plot=Community_Changes1, width=5.5, height=7.1)- ITS2-Mortality
Community_Changes2 <- ggplot (log_Log.data, aes(logC, logD, gropup=Colony, colour=ITS2_14)) +
geom_abline(intercept=0, slope=1, linetype=1, color="black")+
geom_abline(intercept=-1, slope=1, linetype=2, color="black")+
geom_abline(intercept=-2, slope=1, linetype=2, color="darkgray")+
geom_abline(intercept=-3, slope=1, linetype=2, color="gray")+
geom_abline(intercept=-4, slope=1, linetype=2, color="gray")+
geom_abline(intercept=-5, slope=1, linetype=2, color="gray")+
geom_abline(intercept=1, slope=1, linetype=2, color="black")+
geom_abline(intercept=2, slope=1, linetype=2, color="darkgray")+
geom_abline(intercept=3, slope=1, linetype=2, color="gray")+
geom_abline(intercept=4, slope=1, linetype=2, color="gray")+
geom_abline(intercept=5, slope=1, linetype=2, color="gray")+
scale_colour_manual(values=c("blue","purple", "red"),
guide = guide_legend(title.position = "top", nrow = 3),
name = "Initial ITS2",
labels=c("Cd1", "C1b-c", "D1"))+
#scale_colour_manual(values=c("skyblue1", "blue2", "blue4",
# "darkorchid1","darkorchid3",
# "red1", "red3", "red4"))+
ggthe_bw +
theme(legend.position="bottom",
legend.box = "vertical")+
#theme(legend.position="none")+
scale_x_continuous(limits = c(-6,-0.1),
"Cladocopium abundance (log10 C/H))")+
scale_y_continuous(limits = c(-6,-0.1),
"Durusdinium abundance (log10 D/H))") +
geom_point(aes(shape=Mortality),
alpha=0.5, size=3) +
geom_path(size = 0.5 , arrow =
arrow(type = "open", angle = 30,
length = unit(0.1, "inches")))+
annotate("text", x = c(seq(-5.1, -1.1, by=1)), y = -0.9, label = c(
"10,0000D : 1C", "1,000D : 1C", "100D : 1C","10D : 1C", "C : D"), angle = 45)+
annotate("text", x = -0.9, y = c(seq(-4.7, -1.7, by=1)), label = c(
"1D : 10,000C", "1D : 1,000C", "1D : 100C", "1D : 10C"), angle = 45)
Community_Changes2
#ggsave(file="LogLog_labels.svg", plot=Community_Changes2, width=5.5, height=6)- Colony
Community_Changes3 <- ggplot (log_Log.data, aes(logC, logD, gropup=Colony, colour=Colony)) +
geom_abline(intercept=0, slope=1, linetype=1, color="black")+
geom_abline(intercept=-1, slope=1, linetype=2, color="black")+
geom_abline(intercept=-2, slope=1, linetype=2, color="darkgray")+
geom_abline(intercept=-3, slope=1, linetype=2, color="gray")+
geom_abline(intercept=-4, slope=1, linetype=2, color="gray")+
geom_abline(intercept=-5, slope=1, linetype=2, color="gray")+
geom_abline(intercept=1, slope=1, linetype=2, color="black")+
geom_abline(intercept=2, slope=1, linetype=2, color="darkgray")+
geom_abline(intercept=3, slope=1, linetype=2, color="gray")+
geom_abline(intercept=4, slope=1, linetype=2, color="gray")+
geom_abline(intercept=5, slope=1, linetype=2, color="gray")+
#scale_colour_manual(values=c("blue","purple", "red"),
# guide = guide_legend(title.position = "top", nrow = 3),
# name = "Detected symbionts (2014)",
# labels=c("Cladocopium", "Mixed", "Durusdinium"))+
ggthe_bw + theme(legend.position="none")+
scale_x_continuous(limits = c(-6,-0.1), "log10 (C/H cell ratio)")+
scale_y_continuous(limits = c(-6,-0.1), "log10 (D/H cell ratio)") +
geom_point(aes(shape=Yr), alpha=0.5, size=4) +
geom_path(size = 1,arrow =
arrow(type = "open", angle = 30, length = unit(0.1, "inches")))+
annotate("text", x = c(seq(-5.1, -1.1, by=1)), y = -0.9, label = c(
"10,0000D : 1C", "1,000D : 1C", "100D : 1C","10D : 1C", "C : D"), angle = 45)+
annotate("text", x = -0.9, y = c(seq(-4.7, -1.7, by=1)), label = c(
"1D : 10,000C", "1D : 1,000C", "1D : 100C", "1D : 10C"), angle = 45)
Community_Changes3
# ggsave(file="LogLog_Colony.svg", plot=Community_Changes3, width=6, height=6)5c. Changes in the number of colonies dominated by each genus considering mortality
- Fisher test including initially (Aug 2014) 50/50 colonies. Both colonies placed in the D-dominated category to be conservative.
# Only type 1 colonies (n=23)
Input =("
Year C D
2014 9 14
2016 2 18
")
Matriz <- as.matrix(read.table(textConnection(Input),
header=TRUE,
row.names=1))
fisher.test(Matriz,
alternative="two.sided")##
## Fisher's Exact Test for Count Data
##
## data: Matriz
## p-value = 0.03936
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.9354645 61.0268633
## sample estimates:
## odds ratio
## 5.560207# Type 1 and type 3 colonies
Input =("
Year C D
2014 15 14
2016 4 18
")
Matriz <- as.matrix(read.table(textConnection(Input),
header=TRUE,
row.names=1))
fisher.test(Matriz,
alternative="two.sided")##
## Fisher's Exact Test for Count Data
##
## data: Matriz
## p-value = 0.01988
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 1.145585 23.762126
## sample estimates:
## odds ratio
## 4.670456# Spp changes
Input =("
Year I II
2014 23 6
2016 20 2
")
Matriz <- as.matrix(read.table(textConnection(Input),
header=TRUE,
row.names=1))
fisher.test(Matriz,
alternative="two.sided")##
## Fisher's Exact Test for Count Data
##
## data: Matriz
## p-value = 0.4399
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.03477518 2.50826671
## sample estimates:
## odds ratio
## 0.3901317